Análisis Causal de la Mortalidad de Renacuajos de Rana Carrizo

Métodos Analíticos

Introducción

En este proyecto empleamos los datos experimentales de Vonesh & Bolker (2005) [1], quienes en su investigación examinaron las consecuencias de la plasticidad de eclosión inducida por depredadores desde la etapa larval hasta la metamorfosis en la rana de caña de África Oriental, Hyperolius spinigularis [2] realizando un experimento en el que manipulaban el tamaño y la densidad larvaria inicial (imitando los efectos de los depredadores de los huevos). Esperaban que las crías inducidas por depredadores (porque están menos desarrolladas y son más pequeñas) experimentaran mayores tasas de depredación per cápita y un período larvario más largo y, por lo tanto, exhibirían una menor supervivencia a la metamorfosis en presencia de depredadores acuáticos que las larvas más grandes, más desarrolladas y eclosionadas más tarde. Sin embargo, los resultados mostraron que las larvas inducidas por depredadores no solo sobrevivieron a la metamorfosis, sino que también tuvieron tasas de crecimiento más rápidas y alcanzaron tamaños más grandes en la metamorfosis. Esto los motivó a desarrollar un modelo parametrizado a partir de experimentos adicionales para explorar si una combinación de mecanismos, crecimiento compensatorio y depredación específica por densidad y tamaño, podría dar lugar a este patrón. Es por eso que con esta introducción, buscamos replicar y entender su trabajo, utilizando un enfoque bayesiano jerárquico para modelar la mortalidad larval y las respuestas compensatorias posteriores. En este sentido, el modelo jerárquico nos permitirá capturar la heterogeneidad entre los tanques de renacuajos y compartir información entre ellos, lo que resulta en estimaciones más robustas y precisas.

[1] Vonesh, J. R., & Bolker, B. M. (2005). Statistical tools for analyzing larval amphibian survival data. Ecology, 86(1), 172-182.

[2] Se refiere a la capacidad de los embriones de la rana para ajustar su desarrollo y eclosión en respuesta a cambios ambientales, como la presencia de depredadores o el secado de su hábitat

Datos

Los datos provienen de la librería de rethinking.Constan de 48 observaciones que representan los tanques de renacuajos clasificados en pequeños, medianos y grandes, dependiendo de la densidad de renacuajos en cada uno. Además, de información sobre la supervivencia (variable binaria) y de la tasa de supervivencia en cada tanque.

Variable	Descripción
density	Densidad inicial de renacuajos
pred	Factor indicador de presencia de depredadores
size	Tamaño de los renacuajos
surv	Número de renacuajos que sobrevivieron
propsurv	Proporción de supervivencia (surv/density)

En este experimento observamos mucha variación en los datos, y no toda se debe al tratamiento experimental (como la presencia de depredadores). Una gran parte de esa variación proviene de factores no medidos, propios de cada entorno donde viven los renacuajos. Podemos imaginar cada fila del dataset como un “tanque”, es decir, un pequeño ambiente experimental que contiene varios renacuajos. Aunque algunos tanques tengan la misma densidad o condiciones aparentes, hay muchas cosas que no estamos midiendo (como temperatura, luz, microalgas, etc.) que también influyen en la tasa de supervivencia. Esto hace que los tanques funcionen como lo que llamamos un conglomerado o cluster. Dentro de cada tanque observamos múltiples renacuajos, por lo que los datos tienen una estructura agrupada. En otras palabras, tenemos medidas repetidas dentro de grupos que son diferentes entre sí.

Para nuestro análisis nos centraremos en surv como variable de respuesta (binomial) frente a density como total de ensayos.

Es importante mencionar, que si usamos el mismo valor base (intercepto) para todos los tanques pooling, estamos ignorando diferencias importantes entre ellos. Esto puede hacer que no detectemos correctamente el efecto de otras variables como la densidad o el predador. Si por el contrario, usamos un intercepto distinto para cada tanque no pooling pero sin compartir información entre ellos, podríamos caer en lo que se llama una “amnesia estadística”: tratamos a cada tanque como si no tuviéramos nada que aprender de los demás. Pero eso no tiene sentido, porque aunque cada tanque es diferente, los datos de uno pueden ayudarnos a entender mejor a los demás.

Por ello empleamos también un modelo bayesiano jerárquico o multinivel o como lo mencionaremos en este proyecto: modelo de partial pooling con interceptos variables, de este modo, cada tanque tiene su propio parámetro de línea base, y al mismo tiempo estimamos la dispersión entre tanques mediante un prior adaptativo, que aprende de los datos. Con esto buscamor lograr un equilibrio entre asumir que todo es igual (subajuste) y asumir que todo es completamente distinto (sobreajuste).

Nuestros objetivos son:

1. Reproducir y comprender el ejemplo de Statistical Rethinking aplicados a los datos de Reedfrogs.
2. Modelar la mortalidad larval y las respuestas compensatorias posteriores explorando distintos niveles de agrupamiento pooling, no pooling y partial pooling.
3. Evaluar la calidad y complejidad de cada modelo mediante diagnósticos MCMC y criterios de comparación predictiva (WAIC/LOO).
4. Explorar el trade-off underfitting/overfitting mediante simulaciones con distintos tamaños de muestra, ilustrando los beneficios del pooling parcial.
5. Desplegar un análisis causal formal con un DAG que recoja nuestros supuestos de identificación.

Con este enfoque buscamos profundizar en los costes y beneficios de la eclosión temprana inducida por depredadores, y demostrar cómo la regularización adaptativa de los modelos jerárquicos permite inferir efectos individuales de forma más robusta en presencia de datos jerarquizados y dispersos.

DAG

A continuación, se presenta la Gráfica Dirigida Acíclica (DAG) que ilustra las relaciones causales entre las variables de interés. Este DAG se basa en la premisa de que la densidad de renacuajos, el tamaño y la presencia de depredadores influyen en la supervivencia, y que la jerarquía de los tanques también afecta a estas relaciones.

DAG causal: densidad, tamaño, depredadores y jerarquía de tanque

Con

• \(\textrm{T}=\textrm{Tanque}\)

• \(\textrm{D}=\textrm{Densidad inicial}\)

• \(\textrm{G}=\textrm{Tamaño}\)

• \(\textrm{P}=\textrm{Depredadores}\)

• \(\textrm{S}=\textrm{Supervivencia}\)

Aunque podríamos imaginar otras dependencias entre las variables, hay que tener presente que estos datos provienen de un experimento controlado, es decir, como se expone en la introducción el experimento es manipulado y cerrado bajo las condiciones que Vonesh & Bolker establecieron. En un escenario natural, sería razonable investigar vínculos, como el efecto del tamaño de los renacuajos en la densidad poblacional, la influencia de los depredadores sobre esa densidad, o el papel de variables no registradas —por ejemplo, la disponibilidad de alimento u otros recursos— tanto en el tamaño como en la densidad, e incluso factores genéticos que modulen el desarrollo de los renacuajos. Sería muy valioso repetir estas estimaciones en cuerpos de agua naturales, en lugar de en tanques de laboratorio. Por ahora, este experimento nos permite centrarnos en la tasa de supervivencia bajo condiciones estrictamente controladas, estableciendo una base sólida para futuros estudios en la naturaleza.

Modelos

1. Modelo Total Pooling

En este primer modelo totalmente agrupado asumimos que todos los tanques tienen la misma probabilidad de supervivencia. No hay diferencias entre tanques, salvo la variación por la densidad inicial \(D_i\).

\[ S_i \sim \textrm{Binomial}(D_i,p_i) \]

\[ \textrm{logit}(p_i) = \alpha \]

\[ \alpha = \textrm{Normal}(0, 1.5) \] Este modelo ignora la heterogeneidad entre tanques (total pooling) y servirá como línea base para comparar con el modelo jerárquico.

En el enfoque de total pooling, asignamos una única \(\alpha\) a todos los tanques, de modo que el modelo asume idéntica probabilidad de supervivencia para los renacuajos en cada uno de ellos. Es como si tratáramos todos los tanques como réplicas exactas, sin captar ninguna heterogeneidad más allá de la variación provocada por la densidad inicial.

A continuación, mostramos los resultados obtenidos para la estimación de \(\alpha\) en nuestro modelo inicial.

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
alpha	0.84	0.06	0.72	0.84	0.97	1,324.26	1.00
lp__	−685.66	0.70	−687.56	−685.39	−685.17	1,864.17	1.00

Podemos cotejar nuestras predicciones con los datos originales:

Predicciones del modelo 1 (total pooling)

Podemos ver los datos originales (puntos grises) y en las predicciones (puntos rojos) apreciamos que siguen la media de supervivencia de los tanques (línea discontinua azul), lo que evidencia un subajuste del modelo.

2. Modelo No-Pooling

En este modelo no agrupado (no pooling) asignamos un intercepto \(\alpha_i\) distinto a cada tanque, pero no compartimos información entre ellos. Esto significa que cada tanque tiene su propio parámetro de línea base, y no hay aprendizaje entre ellos.

\[ S_i \sim \textrm{Binomial}(D_i,p_i) \]

\[ \textrm{logit}(p_i) = \alpha_i \]

\[ \alpha_i = \textrm{Normal}(0, 1.5)\]

Lo cual se traduce en el siguiente código de Stan:

Podemos ver los resultados de el modelo No-poolong para la estimación de \(\alpha\).

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
alpha[1]	1.71	0.76	0.34	1.66	3.30	7,608.14	1.00
alpha[2]	2.41	0.91	0.80	2.34	4.38	5,328.66	1.00
alpha[3]	0.76	0.63	−0.44	0.74	2.06	7,389.04	1.00
alpha[4]	2.39	0.89	0.79	2.33	4.29	4,871.51	1.00
alpha[5]	1.69	0.76	0.33	1.65	3.27	5,349.81	1.00
alpha[6]	1.72	0.77	0.36	1.67	3.37	4,943.53	1.00
alpha[7]	2.41	0.90	0.84	2.35	4.38	5,784.55	1.00
alpha[8]	1.72	0.76	0.36	1.68	3.35	5,850.42	1.00
alpha[9]	−0.37	0.62	−1.60	−0.37	0.84	6,776.03	1.00
alpha[10]	1.72	0.76	0.36	1.67	3.31	6,189.05	1.00
alpha[11]	0.75	0.62	−0.45	0.73	2.01	6,659.53	1.00
alpha[12]	0.36	0.59	−0.80	0.36	1.53	7,586.50	1.00
alpha[13]	0.76	0.63	−0.42	0.74	1.99	6,438.97	1.00
alpha[14]	0.01	0.58	−1.12	0.01	1.19	6,924.07	1.00
alpha[15]	1.73	0.77	0.32	1.68	3.39	5,612.54	1.00
alpha[16]	1.71	0.77	0.33	1.66	3.38	6,190.12	1.00
alpha[17]	2.55	0.68	1.33	2.50	4.05	5,469.49	1.00
alpha[18]	2.12	0.59	1.06	2.09	3.36	7,128.89	1.00
alpha[19]	1.82	0.54	0.82	1.78	2.96	5,274.61	1.00
alpha[20]	3.10	0.80	1.74	3.04	4.82	4,824.31	1.00
alpha[21]	2.15	0.60	1.06	2.11	3.46	5,935.77	1.00
alpha[22]	2.14	0.62	1.07	2.09	3.45	5,608.99	1.00
alpha[23]	2.13	0.59	1.06	2.09	3.45	6,458.33	1.00
alpha[24]	1.55	0.50	0.67	1.52	2.60	5,688.61	1.00
alpha[25]	−1.11	0.45	−2.05	−1.09	−0.28	5,876.34	1.00
alpha[26]	0.07	0.39	−0.68	0.07	0.85	8,307.26	1.00
alpha[27]	−1.54	0.52	−2.58	−1.51	−0.63	5,461.82	1.00
alpha[28]	−0.55	0.40	−1.33	−0.54	0.21	7,885.78	1.00
alpha[29]	0.08	0.40	−0.69	0.08	0.88	6,664.33	1.00
alpha[30]	1.30	0.47	0.43	1.28	2.29	5,191.47	1.00
alpha[31]	−0.73	0.43	−1.62	−0.72	0.08	6,709.51	1.00
alpha[32]	−0.39	0.42	−1.22	−0.39	0.41	5,782.71	1.00
alpha[33]	2.85	0.66	1.70	2.78	4.30	5,675.78	1.00
alpha[34]	2.46	0.56	1.47	2.43	3.63	5,662.31	1.00
alpha[35]	2.45	0.57	1.44	2.41	3.68	5,372.78	1.00
alpha[36]	1.91	0.48	1.02	1.87	2.95	5,598.11	1.00
alpha[37]	1.91	0.49	1.03	1.88	2.94	5,565.39	1.00
alpha[38]	3.36	0.78	2.05	3.29	5.04	5,814.00	1.00
alpha[39]	2.46	0.58	1.43	2.42	3.73	5,036.40	1.00
alpha[40]	2.16	0.51	1.25	2.13	3.22	6,162.79	1.00
alpha[41]	−1.91	0.48	−2.92	−1.89	−1.03	5,154.94	1.00
alpha[42]	−0.64	0.34	−1.34	−0.63	0.03	5,777.91	1.00
alpha[43]	−0.52	0.34	−1.20	−0.51	0.14	7,900.15	1.00
alpha[44]	−0.40	0.34	−1.06	−0.40	0.25	6,790.40	1.00
alpha[45]	0.51	0.34	−0.16	0.50	1.20	5,771.50	1.00
alpha[46]	−0.63	0.34	−1.32	−0.63	0.02	6,027.82	1.00
alpha[47]	1.92	0.48	1.03	1.88	2.93	6,134.62	1.00
alpha[48]	−0.05	0.33	−0.71	−0.06	0.58	6,124.01	1.00
lp__	−524.22	4.95	−534.85	−523.83	−515.67	1,443.62	1.00

Podemos cotejar nuestras predicciones de este modelo No-pooling con los datos originales:

Predicciones del modelo 2 (no pooling)

Al evaluar las estimaciones del modelo No- pooling (un intercepto \(\alpha_i\) independiente por tanque) obtenemos:

1. Gran variabilidad entre tanques
   • Los \(\alpha_i\) (en escala log-odds) oscilan aproximadamente entre \(-1.9\) y \(+3.4\).
     
   • Transformados a probabilidad, algunos tanques se estiman con supervivencias cercanas al 10–20% y otros al 90–95%.
2. Incertidumbre muy desigual
   • Tanques con pocas observaciones (densidad pequeña) presentan desviaciones estándar de \(\alpha_i\) de 0.6–0.8 y rangos de credibilidad muy amplios (p. ej. \(\alpha_3: [–0.41, +2.06]\).
     
   • Tanques con más renacuajos reducen su incertidumbre a 0.3–0.5 en la desviación estándar de \(\alpha_i\).
3. Sobreajuste
   • Las predicciones del modelo (puntos rojos) siguen casi exactamente los datos observados (puntos grises), incluso en valores extremos.
     
   • No existe “arrastre” hacia un promedio general: cada tanque se ajusta únicamente con su propia información.
4. Problemas en tanques pequeños
   • Con muestras muy pequeñas, pocas muertes o supervivencias cambian drásticamente la estimación de \(\alpha_i\).
     
   • El ancho de los intervalos de credibilidad hace poco útiles esas predicciones para la toma de decisiones.

Por lo anterior, el modelo No-pooling captura fielmente cada dato empírico, pero padece de sobreajuste y de alta incertidumbre en tanques con pocas observaciones. Para obtener estimaciones más estables y evitar extremos sin fundamento, es recomendable utilizar un modelo jerárquico (partial pooling) que comparta información entre tanques.

3. Modelo Partial Pooling

En este tercer modelo parcialmente agrupado (partial pooling) asignamos un intercepto distinto a cada tanque, pero también estimamos la variabilidad entre ellos. Esto nos permite captar la heterogeneidad entre tanques y, al mismo tiempo, compartir información entre ellos. Para estom se agregan dos parámetros: \(\mu\) y \(\sigma\), los cuáles llamaremos hiperparámetros desde este punto.El punto es que, en el modelo anterior, todas las \(\alpha_i\) se distribuían Normal con una media y desviación estándar establecida o fija. Con esta modificación, los \(\alpha_i\) comparten una misma distribución, lo que le permite transmitir información entre tanques al modelo.

\[ S_i \sim \textrm{Binomial}(D_i,p_i) \]

\[ \textrm{logit}(p_i) = \alpha_{T[i]} \]

\[ \alpha_j = \textrm{Normal}(\mu, \sigma) \] \[ \mu \sim \textrm{Normal}(0, 1.5) \] \[ \sigma \sim \textrm{Exponential}(1) \]

Lo cual se traduce en el siguiente código de Stan:

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
alpha_tank[1]	2.15	0.87	0.63	2.08	4.04	4,761.04	1.00
alpha_tank[2]	3.08	1.10	1.19	2.97	5.54	4,249.14	1.00
alpha_tank[3]	1.01	0.67	−0.23	0.99	2.37	5,266.98	1.00
alpha_tank[4]	3.07	1.14	1.16	2.94	5.63	4,224.76	1.00
alpha_tank[5]	2.14	0.89	0.58	2.06	4.04	4,449.29	1.00
alpha_tank[6]	2.14	0.87	0.62	2.07	4.05	4,615.14	1.00
alpha_tank[7]	3.08	1.14	1.17	2.95	5.64	3,843.88	1.00
alpha_tank[8]	2.12	0.86	0.62	2.05	3.99	4,763.29	1.00
alpha_tank[9]	−0.16	0.62	−1.40	−0.15	1.05	6,271.23	1.00
alpha_tank[10]	2.16	0.89	0.60	2.10	4.03	4,241.83	1.00
alpha_tank[11]	1.01	0.66	−0.21	0.99	2.38	5,061.56	1.00
alpha_tank[12]	0.58	0.63	−0.58	0.55	1.89	5,397.86	1.00
alpha_tank[13]	1.01	0.66	−0.19	0.99	2.39	5,147.38	1.00
alpha_tank[14]	0.20	0.61	−1.00	0.19	1.42	5,105.89	1.00
alpha_tank[15]	2.16	0.89	0.63	2.10	4.12	4,749.16	1.00
alpha_tank[16]	2.16	0.87	0.62	2.10	4.09	4,319.13	1.00
alpha_tank[17]	2.89	0.78	1.57	2.82	4.63	4,468.31	1.00
alpha_tank[18]	2.41	0.67	1.24	2.36	3.83	4,593.79	1.00
alpha_tank[19]	2.01	0.60	0.94	1.98	3.33	5,093.56	1.00
alpha_tank[20]	3.65	1.02	1.95	3.54	5.99	4,236.48	1.00
alpha_tank[21]	2.40	0.66	1.24	2.34	3.86	4,093.93	1.00
alpha_tank[22]	2.40	0.67	1.22	2.35	3.88	4,798.60	1.00
alpha_tank[23]	2.38	0.66	1.21	2.34	3.86	4,266.05	1.00
alpha_tank[24]	1.71	0.54	0.75	1.68	2.85	4,827.15	1.00
alpha_tank[25]	−0.99	0.45	−1.91	−0.98	−0.12	5,568.18	1.00
alpha_tank[26]	0.17	0.40	−0.61	0.16	0.94	5,482.17	1.00
alpha_tank[27]	−1.42	0.48	−2.42	−1.39	−0.56	5,050.46	1.00
alpha_tank[28]	−0.46	0.40	−1.24	−0.46	0.30	5,070.96	1.00
alpha_tank[29]	0.16	0.39	−0.60	0.16	0.92	5,649.98	1.00
alpha_tank[30]	1.44	0.48	0.56	1.42	2.43	5,226.81	1.00
alpha_tank[31]	−0.63	0.42	−1.49	−0.63	0.17	5,171.55	1.00
alpha_tank[32]	−0.31	0.39	−1.12	−0.31	0.43	5,537.69	1.00
alpha_tank[33]	3.19	0.77	1.88	3.14	4.81	4,027.19	1.00
alpha_tank[34]	2.71	0.63	1.61	2.67	4.05	4,766.21	1.00
alpha_tank[35]	2.72	0.64	1.59	2.66	4.11	3,616.00	1.00
alpha_tank[36]	2.07	0.52	1.17	2.04	3.15	4,635.30	1.00
alpha_tank[37]	2.05	0.52	1.10	2.03	3.15	5,464.88	1.00
alpha_tank[38]	3.91	0.98	2.31	3.80	6.12	3,915.94	1.00
alpha_tank[39]	2.71	0.66	1.57	2.67	4.13	4,785.29	1.00
alpha_tank[40]	2.35	0.58	1.31	2.31	3.61	4,926.42	1.00
alpha_tank[41]	−1.81	0.47	−2.82	−1.78	−0.94	5,490.64	1.00
alpha_tank[42]	−0.58	0.34	−1.26	−0.58	0.05	5,100.21	1.00
alpha_tank[43]	−0.45	0.34	−1.14	−0.45	0.21	4,885.82	1.00
alpha_tank[44]	−0.33	0.34	−1.01	−0.33	0.32	4,534.55	1.00
alpha_tank[45]	0.58	0.34	−0.07	0.58	1.24	6,064.61	1.00
alpha_tank[46]	−0.57	0.36	−1.30	−0.57	0.13	5,006.21	1.00
alpha_tank[47]	2.07	0.52	1.14	2.05	3.19	4,712.36	1.00
alpha_tank[48]	0.01	0.34	−0.66	0.01	0.68	6,330.02	1.00
lp__	−532.98	5.54	−545.05	−532.57	−523.44	1,487.68	1.00

Y para los hiperparámetros \(\mu\) y \(\sigma\):

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
mu_alpha	1.35	0.26	0.86	1.35	1.85	3,749.40	1.00
sigma_alpha	1.62	0.22	1.24	1.60	2.08	2,568.13	1.00

Podemos nuevamente, ver cómo se comportan nuestras predicciones contra los valores observados:

Predicciones del modelo 3 (partial pooling)

Esta estructura jerárquica permite al modelo “aprender entre tanques”: los tanques con poca información se benefician del conocimiento global, mientras que los tanques con muchos datos retienen sus características particulares. Esto se conoce como regularización adaptativa o “shrinkage”.

Las estimaciones de los hiperparámetros: - La media global \(\mu_\alpha\) ≈ 1.35, lo que se traduce en una tasa de supervivencia promedio de aproximadamente 79% (vía función logística), es decir por arriba del promedio. - La desviación estándar \(\sigma_\alpha\) ≈ 1.62, lo cual indica una heterogeneidad considerable entre tanques.

Los interceptos de tanques con pocos renacuajos (menos información) se encogen más hacia la media global, generando estimaciones más conservadoras.En contraste, los tanques con mayor número de observaciones muestran menos shrinkage, permitiendo conservar la variación real.

El modelo Parcial pooling corrige el subajuste del modelo Total pooling y el sobreajuste del modelo No-pooling, logrando un balance óptimo.Como podemos ver en la última gráfica, las predicciones siguen de cerca los datos observados, pero sin reaccionar de forma extrema a variaciones pequeñas o poco informativas. También podemos observar en la tabla de resultados que los valores de \(\hat{R}\) ≈ 1.00 y los tamaños efectivos (n_eff) son altos, lo cual indica que las cadenas MCMC convergieron adecuadamente y los resultados son confiables.

Con esto, si podemos concluir que el modelo jerárquico permite identificar patrones globales sin ignorar las diferencias individuales, siendo especialmente útil cuando algunos grupos (tanques) tienen poca información, con esto confirmas que existe variabilidad realen las tasas de supervivencia entre tanques.

Relación entre log-odds y probabilidad

Para comprender mejor el efecto del parámetro α en los modelos logísticos, se muestra a continuación la transformación inv_logit(α):

Esto ilustra cómo pequeñas diferencias en log-odds (α) se traducen en cambios más o menos pronunciados en la probabilidad, dependiendo de la región de la curva sigmoide.

Comparativos

Podemos graficar las estimaciones del modelo No-pooling y el Parcial pooling para ver el comportamiento con mayor detenimiento.

Incluyendo intervalos de credibilidad para ver su comportamiento:

Ahora que tenemos los tres modelos ajustados, es momento de compararlos. Para ello, utilizaremos: - Diagnósticos MCMC: Para confirmar la buena convergencia de los modelos. - Visualización de predicciones: Para contrastar la calidad de las predicciones por tanque.

Diagnósticos MCMC

Una buena práctica es comparar las trazas y distribuciones posteriores de parámetros clave para verificar convergencia y comportamiento estable.

Estos gráficos nos permiten confirmar que las cadenas exploran bien el espacio posterior, sin signos de problemas como falta de mezcla o divergencias.

Predicciones por tanque

Unificamos en una sola gráfica las predicciones de los tres modelos para compararlas directamente.

Este gráfico muestra cómo cada modelo aproxima los datos observados. Veremos que el Modelo parcial poolong logra un buen balance, sin subestimar ni sobreajustar.

Modelo con información completa

En los tres modelos anteriores no incorporamos todas las variables, y como lo que deseamos es estudiar el fenómeno, incluir todas las variables, es decir, agregando la variable depredador P y la variable tamaño G, dado nuestro DAG, puede mejorar la estimación de los efectos por variable. Podemos extender entonces la idea del Modelo Jerárquico e incorporar estas variables. De tal modo que:

\[ S_i \sim \textrm{Binomial}(D_i, p_i) \]

\[ \textrm{logit}(p_i) = \alpha_{\textrm{Tanque}[i]} + \beta_p *\textrm{pred} + \beta_s *\textrm{size}\]

\[ \alpha_j \sim \textrm{Normal}(\mu, \sigma) \]

\[ \mu \sim \textrm{Normal}(0, 1.5) \]

\[ \sigma \sim \textrm{Exponential}(1) \]

\[ \beta_p \sim \textrm{Normal}(-0.5,1) \]

\[ \beta_s \sim \textrm{Normal}(0,1) \]

Por lo que el codigo en stan:

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
alpha_tank[1]	2.74	0.69	1.38	2.73	4.15	2,585.31	1.00
alpha_tank[2]	3.16	0.71	1.81	3.14	4.62	2,445.93	1.00
alpha_tank[3]	2.02	0.64	0.71	2.02	3.30	1,557.38	1.00
alpha_tank[4]	3.15	0.72	1.79	3.13	4.65	2,262.92	1.00
alpha_tank[5]	2.61	0.68	1.37	2.61	4.00	2,758.04	1.00
alpha_tank[6]	2.61	0.66	1.33	2.60	3.96	3,176.61	1.00
alpha_tank[7]	3.06	0.72	1.74	3.03	4.53	3,185.82	1.00
alpha_tank[8]	2.61	0.69	1.31	2.59	4.00	3,820.07	1.00
alpha_tank[9]	2.54	0.61	1.28	2.55	3.69	952.45	1.00
alpha_tank[10]	3.80	0.64	2.59	3.78	5.11	1,304.93	1.01
alpha_tank[11]	3.27	0.62	2.06	3.26	4.47	1,126.72	1.00
alpha_tank[12]	3.02	0.61	1.84	3.00	4.21	1,011.67	1.00
alpha_tank[13]	3.04	0.58	1.91	3.04	4.20	1,714.41	1.00
alpha_tank[14]	2.54	0.57	1.42	2.53	3.66	1,539.77	1.00
alpha_tank[15]	3.60	0.61	2.47	3.58	4.84	1,919.14	1.00
alpha_tank[16]	3.60	0.63	2.42	3.56	4.96	2,222.43	1.00
alpha_tank[17]	3.15	0.61	2.01	3.13	4.41	2,282.52	1.00
alpha_tank[18]	2.85	0.57	1.72	2.84	3.99	1,935.06	1.00
alpha_tank[19]	2.59	0.54	1.56	2.60	3.69	1,457.11	1.01
alpha_tank[20]	3.48	0.66	2.30	3.44	4.87	2,440.63	1.00
alpha_tank[21]	2.65	0.56	1.61	2.63	3.78	3,979.88	1.00
alpha_tank[22]	2.64	0.56	1.60	2.62	3.77	3,741.02	1.00
alpha_tank[23]	2.65	0.57	1.59	2.65	3.83	3,521.71	1.00
alpha_tank[24]	2.11	0.52	1.12	2.10	3.12	3,533.69	1.00
alpha_tank[25]	1.92	0.55	0.76	1.95	2.96	739.07	1.01
alpha_tank[26]	2.86	0.52	1.82	2.87	3.84	640.67	1.01
alpha_tank[27]	1.61	0.58	0.42	1.63	2.71	746.12	1.01
alpha_tank[28]	2.35	0.54	1.26	2.35	3.37	671.43	1.01
alpha_tank[29]	2.54	0.46	1.64	2.55	3.44	1,131.69	1.00
alpha_tank[30]	3.51	0.48	2.62	3.50	4.47	1,392.06	1.00
alpha_tank[31]	1.90	0.47	0.95	1.91	2.78	1,134.74	1.00
alpha_tank[32]	2.16	0.46	1.26	2.17	3.04	1,144.47	1.00
alpha_tank[33]	3.33	0.58	2.22	3.31	4.53	2,150.16	1.00
alpha_tank[34]	3.05	0.56	2.01	3.03	4.23	1,742.16	1.00
alpha_tank[35]	3.04	0.55	2.04	3.03	4.18	1,585.32	1.01
alpha_tank[36]	2.60	0.52	1.64	2.60	3.64	1,380.73	1.00
alpha_tank[37]	2.33	0.48	1.46	2.31	3.31	3,860.05	1.00
alpha_tank[38]	3.47	0.64	2.34	3.43	4.80	3,177.58	1.00
alpha_tank[39]	2.83	0.56	1.82	2.80	3.99	4,385.16	1.00
alpha_tank[40]	2.56	0.51	1.60	2.54	3.63	4,245.97	1.00
alpha_tank[41]	1.29	0.58	0.09	1.32	2.35	710.92	1.01
alpha_tank[42]	2.26	0.51	1.23	2.28	3.21	604.59	1.01
alpha_tank[43]	2.36	0.51	1.33	2.37	3.34	605.40	1.01
alpha_tank[44]	2.47	0.51	1.46	2.48	3.43	578.77	1.01
alpha_tank[45]	2.91	0.43	2.04	2.91	3.73	1,039.82	1.00
alpha_tank[46]	1.93	0.43	1.06	1.93	2.77	968.73	1.00
alpha_tank[47]	3.99	0.49	3.05	3.98	4.99	1,490.47	1.00
alpha_tank[48]	2.42	0.43	1.55	2.42	3.25	967.64	1.00
lp__	−500.70	6.71	−514.41	−500.46	−488.42	937.50	1.00

Y para los hiperparámetros \(\mu\), \(\sigma\) y las \(\beta\)´s :

parameter	mean	sd	2.5%	50%	97.5%	n_eff	Rhat
mu_alpha	2.72	0.28	2.18	2.72	3.26	462.58	1.01
sigma_alpha	0.78	0.15	0.53	0.77	1.10	1,271.72	1.00
bp	−2.42	0.30	−2.99	−2.42	−1.81	601.71	1.01
bs	−0.41	0.28	−0.96	−0.42	0.17	622.89	1.01

Como en los modelos anteriores, podemos ver cómo se comportan nuestras predicciones de este modelo completo contra los valores observados:

Si comparamos este modelo completo contra el Parcial pooling que era hasta este momento el mejor modelo, podemos visualizar esta comparación:

Pareciera que al incorporar más información el modelo no lo captura ya que se alejan un poco más de lo observado, si observamos la parte de mayor tasa de supervivencia en ambas gráficas, los datos de los modelos se alejan más de los observados, lo que podría indicar que el modelo no está capturando la información de los depredadores y el tamaño de los renacuajos.

Conclusiones

En este proyecto aplicamos modelos jerárquicos bayesianos para analizar datos experimentales sobre la supervivencia de renacuajos, recolectados por Vonesh & Bolker en un entorno controlado con variaciones en densidad poblacional y presencia de depredadores. A través de un enfoque comparativo, evaluamos tres estrategias de modelación —total pooling, no pooling y partial pooling— y exploramos sus implicaciones estadísticas y ecológicas.

Nuestros principales hallazgos son:

Los modelos Total pooling subestiman la variabilidad

El enfoque de total pooling ignora la estructura jerárquica de los datos, asumiendo una única probabilidad de supervivencia para todos los tanques. Esto conduce a un subajuste (underfitting), perdiendo de vista diferencias importantes debidas a factores no observados, como la variación ambiental o genética entre tanques.

Los modelos No-pooling sobreajustan

El modelo sin agrupamiento (no pooling) asigna parámetros independientes a cada tanque, sin compartir información entre ellos. Esto resulta en sobreajuste (overfitting), especialmente en tanques con pocos datos, donde incluso pequeñas fluctuaciones generan estimaciones inestables.

El modelo jerárquico o Parcial pooling logra un equilibrio

El modelo parcialmente agrupado (partial pooling) combina la flexibilidad de asignar parámetros específicos a cada grupo con la estabilidad de una distribución común que los regula. Esta estructura: - Captura la heterogeneidad real entre tanques. - Reduce la incertidumbre en grupos pequeños mediante shrinkage adaptativo. - Mejora la capacidad predictiva y generalización del modelo.

A nivel ecológico, encontramos que:

• La supervivencia promedio fue mayor en renacuajos pequeños, un resultado que, aunque contraintuitivo, coincide con las conclusiones del estudio original.
• La densidad inicial afecta negativamente la supervivencia, probablemente por competencia intraespecífica o mayor exposición a depredadores.
• La presencia de depredadores incrementa la variabilidad de las tasas de supervivencia, posiblemente por efectos de habituación o respuestas defensivas plásticas.
• Las distribuciones iniciales de supervivencia pueden actuar como regularizadores informativos, permitiendo una inferencia más robusta.